function [mask newimg] = KidneyMask(vol,sigma, pixfrom,pixto,factor)
%% Filter kidney image and show mask.  Calculate average along z-axis
% in region not containing sample.
%  Input args:
%   vol:     3D volume
%   sigma:   standard deviation of Gaussian filter
%   pixfrom: index starting the sample averaging (must be in range of z-axis elements)
%   pixto:   index ending the averaging sample
%  Output: 
%   mask:    2D image showing the filtered background average of the image 
%   newimg:  new volume that has been filtered 
%
%   - Michael Eager,  (michael.eager@monash.edu)
%   - (c) 2012,  Monash Biomedical Imaging,  Monash University, Australia


%corr = inhomogenityCorrection(vol,sigma);
len = size(vol);

% Multiplcation factor to mask
if nargin < 5
factor=5;
end

% Calculate background along main axis
bg =  squeeze(mean(vol(pixfrom:pixto,:,:),1));

% Smooth background image
f = fspecial('gaussian', 10, sigma);
mask = imfilter(bg,f,'replicate');

% Normalise the peak value to 1
mask = mask./max(mask(:));

for ii = 1:len(1)
        newimg(:,:,ii) = squeeze(vol(ii,:,:))./(factor.*mask);
end

% reorganise volume and normalise it
newimg = permute(newimg,[3 1 2]);
newimg= NormaliseImage(newimg);
